Additional Figures

Accompanying ‘CONVOCALS_downstream’

General Setup

General Setup

Create a new analysis directories.

- general directory

- stains directory

- for plots

- for output of summary results

- for baseline tables

- for genetic analyses

- for Cox regression results

source("scripts/functions.R")

source("scripts/pack02.packages.R")

* General packages...

* Genomic packages...

Today = format(as.Date(as.POSIXlt(Sys.time())), "%Y%m%d")
Today.Report = format(as.Date(as.POSIXlt(Sys.time())), "%A, %B %d, %Y")
source("scripts/colors.R")

Background

This notebook contains additional figures of the project.

Loading data

load(paste0(PROJECT_loc, "/",Today,".",PROJECTNAME,".bulkRNAseq.main_analysis.RData"))
# load(paste0(PROJECT_loc, "/20241017.",PROJECTNAME,".bulkRNAseq.main_analysis.RData"))

Fix some variables

We need to get the ‘conventional unit’ versions of cholesterols.

AERNASE.clin <- merge(AERNASE.clin.targets, 
                            subset(AEDB.CEA, select = c("STUDY_NUMBER", 
                                                        "risk614", 
                                                        "LDL_finalCU", "HDL_finalCU", "TC_finalCU", "TG_finalCU")), 
                            by.x = "study_number", by.y = "STUDY_NUMBER", sort = TRUE, all.x = TRUE)

Additional figures

Age and sex

We want to create per-age-group figures median ± interquartile range.

• Box and Whisker plot for CONVOCALS_downstream plaque levels by sex.
• Box and Whisker plot for CONVOCALS_downstream plaque levels by (sex and) age group (<55, 55-64, 65-74, 75-84, 85+).

# ?ggpubr::ggboxplot()
compare_means(HMOX1 ~ Gender,  data = AERNASE.clin, method = "kruskal.test")

ggpubr::ggboxplot(AERNASE.clin, 
                  x = c("Gender"),
                  y = "HMOX1", 
                  xlab = "gender",
                  ylab = "HMOX1 (normalized expression)",
                  color = "Gender",
                  palette = c("#D5267B", "#1290D9"),
                  add = "jitter") +
  stat_compare_means(aes(group = Gender), label = "p.format", method = "kruskal.test")
ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.Gender.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

library(dplyr)

AERNASE.clin <- AERNASE.clin %>% dplyr::mutate(AgeGroup = factor(case_when(Age < 55 ~ "<55",
                                                     Age >= 55  & Age <= 64 ~ "55-64",
                                                     Age >= 65  & Age <= 74 ~ "65-74",
                                                     Age >= 75  & Age <= 84 ~ "75-84",
                                                     Age >= 85 ~ "85+"))) 

AERNASE.clin <- AERNASE.clin %>% dplyr::mutate(AgeGroupSex = factor(case_when(Age < 55 & Gender == "male" ~ "<55 males" ,
                                                        Age >= 55  & Age <= 64 & Gender == "male"~ "55-64 males",
                                                        Age >= 65  & Age <= 74 & Gender == "male"~ "65-74 males",
                                                        Age >= 75  & Age <= 84 & Gender == "male"~ "75-84 males",
                                                        Age >= 85 & Gender == "male"~ "85+ males",
                                                        Age < 55 & Gender == "female" ~ "<55 females" ,
                                                        Age >= 55  & Age <= 64 & Gender == "female"~ "55-64 females ",
                                                        Age >= 65  & Age <= 74 & Gender == "female"~ "65-74 females",
                                                        Age >= 75  & Age <= 84 & Gender == "female"~ "75-84 females",
                                                        Age >= 85 & Gender == "female"~ "85+ females")))

table(AERNASE.clin$AgeGroup, AERNASE.clin$Gender)

       
        female male
  <55       16   38
  55-64     82  195
  65-74    114  318
  75-84     84  213
  85+       12   20

table(AERNASE.clin$AgeGroupSex)

   <55 females      <55 males 55-64 females     55-64 males  65-74 females    65-74 males  75-84 females    75-84 males    85+ females      85+ males 
            16             38             82            195            114            318             84            213             12             20 

Now we can draw some graphs of plaque CONVOCALS_downstream levels per sex and age group as median ± interquartile range.

# ?ggpubr::ggboxplot()
compare_means(HMOX1 ~ AgeGroup,  data = AERNASE.clin, method = "kruskal.test")

ggpubr::ggboxplot(AERNASE.clin, 
                  x = c("AgeGroup"),
                  y = "HMOX1", 
                  xlab = "Age groups (years)",
                  ylab = "HMOX1 (normalized expression)",
                  color = "AgeGroup",
                  palette = "npg",
                  # add = "median_iqr")
                  add = c("median_iqr", "jitter")) +
  stat_compare_means(aes(group = AgeGroup), label = "p.format", method = "kruskal.test")
ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.AgeGroup.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ AgeGroup, group.by = "Gender", data = AERNASE.clin, method = "kruskal.test")

ggpubr::ggboxplot(AERNASE.clin, 
                  x = c("AgeGroup"),
                  y = "HMOX1", 
                  xlab = "Age groups (years) per gender",
                  ylab = "HMOX1 (normalized expression",
                  color = "Gender",
                  palette = c("#D5267B", "#1290D9"),
                  # add = "median_iqr")
                  add = c("median_iqr", "jitter")) +
  stat_compare_means(aes(group = Gender), label = "p.format", method = "kruskal.test")
ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.AgeGroup_perGender.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

Hypertension & blood pressure

We want to create figures of CONVOCALS_downstream levels stratified by hypertension/blood pressure, and use of anti-hypertensive drugs.

• Box and Whisker plot for CONVOCALS_downstream plaque levels by hypertension group (no, yes)
• Box and Whisker plot for CONVOCALS_downstream plaque levels by systolic blood pressure group (<120, 120-139, 140-159,160+)

library(dplyr)

AERNASE.clin <- AERNASE.clin %>% mutate(SBPGroup = factor(case_when(systolic < 120 ~ "<120",
                                                     systolic >= 120  & systolic <= 139 ~ "120-139",
                                                     systolic >= 140  & systolic <= 159 ~ "140-159",
                                                     systolic >= 160 ~ "160+"))) 

table(AERNASE.clin$SBPGroup, AERNASE.clin$Gender)

         
          female male
  <120        18   57
  120-139     63  138
  140-159     79  211
  160+       116  260

Now we can draw some graphs of plaque CONVOCALS_downstream levels per sex and hypertension/blood pressure group as median ± interquartile range.

HMOX1

detach("package:EnsDb.Hsapiens.v86", unload = TRUE)
detach("package:ensembldb", unload = TRUE)

compare_means(HMOX1 ~ SBPGroup, data = AERNASE.clin %>% filter(!is.na(SBPGroup)), method = "kruskal.test")

filter: removed 150 rows (14%), 942 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(SBPGroup)), 
                  x = c("SBPGroup"),
                  y = "HMOX1", 
                  xlab = "Systolic blood pressure (mmHg)",
                  ylab = "HMOX1 (normalized expression)",
                  color = "SBPGroup",
                  palette = "npg",
                  add = "jitter") +
  stat_compare_means(aes(group = SBPGroup), label = "p.format", method = "kruskal.test")

filter: removed 150 rows (14%), 942 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.SBPGroup.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ Hypertension.selfreport, data = AERNASE.clin %>% filter(!is.na(Hypertension.selfreport)), method = "kruskal.test")

filter: removed 30 rows (3%), 1,062 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(Hypertension.selfreport)), 
                  x = c("Hypertension.selfreport"),
                  y = "HMOX1", 
                  xlab = "Self-reported hypertension",
                  ylab = "HMOX1 (normalized expression)",
                  color = "Hypertension.selfreport",
                  palette = "npg",
                  add = "jitter") +
  stat_compare_means(aes(group = Hypertension.selfreport), label = "p.format", method = "kruskal.test")

filter: removed 30 rows (3%), 1,062 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.Hypertension.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ Hypertension.drugs, data = AERNASE.clin %>% filter(!is.na(Hypertension.drugs)), method = "kruskal.test")

filter: removed one row (<1%), 1,091 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(Hypertension.drugs)), 
                  x = c("Hypertension.drugs"),
                  y = "HMOX1", 
                  xlab = "Hypertension medication use",
                  ylab = "HMOX1 (normalized expression)",
                  color = "Hypertension.drugs",
                  palette = "npg",
                  add = "jitter") +
  stat_compare_means(aes(group = Hypertension.drugs), label = "p.format", method = "kruskal.test")

filter: removed one row (<1%), 1,091 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.HypertensionDrugs.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ SBPGroup, group.by = "Gender", data = AERNASE.clin %>% filter(!is.na(SBPGroup)), method = "kruskal.test")

filter: removed 150 rows (14%), 942 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(SBPGroup)), 
                  x = c("SBPGroup"),
                  y = "HMOX1", 
                  xlab = "Systolic blood pressure (mmHg) per gender",
                  ylab = "HMOX1 (normalized expression)",
                  color = "Gender",
                  palette = c("#D5267B", "#1290D9"),
                  add = "jitter") +
  stat_compare_means(aes(group = Gender), label = "p.format", method = "kruskal.test")

filter: removed 150 rows (14%), 942 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.SBPGroup_byGender.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ Hypertension.selfreport, group.by = "Gender", data = AERNASE.clin %>% filter(!is.na(Hypertension.selfreport)), method = "kruskal.test")

filter: removed 30 rows (3%), 1,062 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(Hypertension.selfreport)), 
                  x = c("Hypertension.selfreport"),
                  y = "HMOX1", 
                  xlab = "Self-reported hypertension per gender",
                  ylab = "HMOX1 (normalized expression)",
                  color = "Gender",
                  palette = c("#D5267B", "#1290D9"),
                  add = "jitter") +
  stat_compare_means(aes(group = Gender), label = "p.format", method = "kruskal.test")

filter: removed 30 rows (3%), 1,062 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.Hypertension_byGender.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ Hypertension.drugs, group.by = "Gender", data = AERNASE.clin %>% filter(!is.na(Hypertension.drugs)), method = "kruskal.test")

filter: removed one row (<1%), 1,091 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(Hypertension.drugs)), 
                  x = c("Hypertension.drugs"),
                  y = "HMOX1", 
                  xlab = "Hypertension medication use per gender",
                  ylab = "HMOX1 (normalized expression)",
                  color = "Gender",
                  palette = c("#D5267B", "#1290D9"),
                  add = "jitter") +
  stat_compare_means(aes(group = Gender), label = "p.format", method = "kruskal.test")

filter: removed one row (<1%), 1,091 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.Hypertension.drugs_byGender.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ SBPGroup, group.by = "Hypertension.drugs", data = AERNASE.clin %>% filter(!is.na(SBPGroup) & !is.na(Hypertension.drugs)), method = "kruskal.test")

filter: removed 151 rows (14%), 941 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(SBPGroup) & !is.na(Hypertension.drugs)), 
                  x = c("SBPGroup"),
                  y = "HMOX1", 
                  xlab = "Systolic blood pressure (mmHg) by medication use",
                  ylab = "HMOX1 (normalized expression)",
                  color = "Hypertension.drugs",
                  palette = c("#49A01D", "#1290D9"),
                  add = "jitter") +
  stat_compare_means(aes(group = Hypertension.drugs), label = "p.format", method = "kruskal.test")

filter: removed 151 rows (14%), 941 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.SBPGroup_byHypertensionDrugs.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ Hypertension.selfreport, group.by = "Hypertension.drugs", data = AERNASE.clin %>% filter(!is.na(Hypertension.selfreport) & !is.na(Hypertension.drugs)), method = "kruskal.test")

filter: removed 31 rows (3%), 1,061 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(Hypertension.selfreport) & !is.na(Hypertension.drugs)), 
                  x = c("Hypertension.selfreport"),
                  y = "HMOX1", 
                  xlab = "Self-reported hypertension per medication use",
                  ylab = "HMOX1 (normalized expression)",
                  color = "Hypertension.drugs",
                  palette = c("#49A01D", "#1290D9"),
                  add = "jitter") +
  stat_compare_means(aes(group = Hypertension.drugs), label = "p.format", method = "kruskal.test")

filter: removed 31 rows (3%), 1,061 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.Hypertension.selfreport_byHypertensionDrugs.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

Hypercholesterolemia & LDL levels

We want to create figures of CONVOCALS_downstream levels stratified by hypercholesterolemia/LDL-levels, and use of lipid-lowering drugs.

• Box and Whisker plot for CONVOCALS_downstream plaque levels by hypercholesterolemia (risk614) group (no, yes)
• Box and Whisker plot for CONVOCALS_downstream plaque levels by lipid-lowering drugs group (no, yes)
• Box and Whisker plot for CONVOCALS_downstream plaque levels by LDL-levels (mmol/L) group (<100, 100-129, 130-159, 160-189, 190+)

library(dplyr)

AERNASE.clin <- AERNASE.clin %>% mutate(LDLGroup = factor(case_when(LDL_finalCU < 100 ~ "<100",
                                                     LDL_finalCU >= 100  & LDL_finalCU <= 129 ~ "100-129",
                                                     LDL_finalCU >= 130  & LDL_finalCU <= 159 ~ "130-159",
                                                     LDL_finalCU >= 160  & LDL_finalCU <= 189 ~ "160-189",
                                                     LDL_finalCU >= 190 ~ "190+"))) 


table(AERNASE.clin$LDLGroup, AERNASE.clin$Gender)

         
          female male
  <100        83  219
  100-129     49  120
  130-159     33   68
  160-189     16   24
  190+         9   14

require(sjlabelled)

AERNASE.clin$risk614 <- to_factor(AERNASE.clin$risk614)

# Fix plaquephenotypes
attach(AERNASE.clin)
AERNASE.clin[,"Hypercholesterolemia"] <- NA
AERNASE.clin$Hypercholesterolemia[risk614 == "missing value"] <- NA
AERNASE.clin$Hypercholesterolemia[risk614 == -999] <- NA
AERNASE.clin$Hypercholesterolemia[risk614 == 0] <- "no"
AERNASE.clin$Hypercholesterolemia[risk614 == 1] <- "yes"
detach(AERNASE.clin)

table(AERNASE.clin$risk614, AERNASE.clin$Hypercholesterolemia)

   
     no yes
  0 306   0
  1   0 697

# AEDB.temp <- subset(AEDB,  select = c("STUDY_NUMBER", "UPID", "Age", "Gender", "Hospital", "Artery_summary", "risk614", "Hypercholesterolemia"))
# require(labelled)
# AEDB.temp$Gender <- to_factor(AEDB.temp$Gender)
# AEDB.temp$Hospital <- to_factor(AEDB.temp$Hospital)
# AEDB.temp$Artery_summary <- to_factor(AEDB.temp$Artery_summary)
# 
# DT::datatable(AEDB.temp[1:10,], caption = "Excerpt of the whole AEDB.", rownames = FALSE)
# 
# rm(AEDB.temp)

Now we can draw some graphs of plaque CONVOCALS_downstream levels per sex and hypercholesterolemia/LDL-levels group, as well as stratified by lipid-lowering drugs users as median ± interquartile range.

HMOX1

compare_means(HMOX1 ~ LDLGroup, data = AERNASE.clin %>% filter(!is.na(LDLGroup)), method = "kruskal.test")

filter: removed 457 rows (42%), 635 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(LDLGroup)), 
                  x = c("LDLGroup"),
                  y = "HMOX1", 
                  xlab = "LDL (mg/dL) per gender",
                  ylab = "HMOX1 (normalized expression))",
                  color = "LDLGroup",
                  palette = "npg",
                  add = "jitter") +
  stat_compare_means(label = "p.format", method = "kruskal.test")

filter: removed 457 rows (42%), 635 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.LDLGroups.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ LDLGroup, group.by = "Gender", data = AERNASE.clin %>% filter(!is.na(LDLGroup)), method = "kruskal.test")

filter: removed 457 rows (42%), 635 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(LDLGroup)), 
                  x = c("LDLGroup"),
                  y = "HMOX1", 
                  xlab = "LDL (mg/dL) per gender",
                  ylab = "HMOX1 (normalized expression))",
                  color = "Gender",
                  palette = c("#D5267B", "#1290D9"),
                  add = "jitter") +
  stat_compare_means(aes(group = Gender), label = "p.format", method = "kruskal.test")

filter: removed 457 rows (42%), 635 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.LDLGroups_byGender.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ Hypercholesterolemia, data = AERNASE.clin %>% filter(!is.na(Hypercholesterolemia)), method = "kruskal.test")

filter: removed 89 rows (8%), 1,003 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(Hypercholesterolemia)), 
                  x = c("Hypercholesterolemia"),
                  y = "HMOX1", 
                  xlab = "Diagnosed hypercholesterolemia",
                  ylab = "HMOX1 (normalized expression))",
                  color = "Hypercholesterolemia",
                  palette = "npg",
                  add = "jitter") +
  stat_compare_means(label = "p.format", method = "kruskal.test")

filter: removed 89 rows (8%), 1,003 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.Hypercholesterolemia.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ Hypercholesterolemia, group.by = "Gender", data = AERNASE.clin %>% filter(!is.na(Hypercholesterolemia)), method = "kruskal.test")

filter: removed 89 rows (8%), 1,003 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(Hypercholesterolemia)), 
                  x = c("Hypercholesterolemia"),
                  y = "HMOX1", 
                  xlab = "Diagnosed hypercholesterolemia per gender",
                  ylab = "HMOX1 (normalized expression))",
                  color = "Gender",
                  palette = c("#D5267B", "#1290D9"),
                  add = "jitter") +
  stat_compare_means(aes(group = Gender), label = "p.format", method = "kruskal.test")

filter: removed 89 rows (8%), 1,003 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.Hypercholesterolemia_byGender.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ Med.Statin.LLD, data = AERNASE.clin %>% filter(!is.na(Med.Statin.LLD)), method = "kruskal.test")

filter: removed one row (<1%), 1,091 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(Med.Statin.LLD)), 
                  x = c("Med.Statin.LLD"),
                  y = "HMOX1", 
                  xlab = "Lipid-lowering drug use",
                  ylab = "HMOX1 (normalized expression))",
                  color = "Med.Statin.LLD",
                  palette = "npg",
                  add = "jitter") +
  stat_compare_means(label = "p.format", method = "kruskal.test")

filter: removed one row (<1%), 1,091 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.Med.Statin.LLD.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ Med.Statin.LLD, group.by = "Gender", data = AERNASE.clin %>% filter(!is.na(Med.Statin.LLD)), method = "kruskal.test")

filter: removed one row (<1%), 1,091 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(Med.Statin.LLD)), 
                  x = c("Med.Statin.LLD"),
                  y = "HMOX1", 
                  xlab = "Lipid-lowering drug use per gender",
                  ylab = "HMOX1 (normalized expression))",
                  color = "Gender",
                  palette = c("#D5267B", "#1290D9"),
                  add = "jitter") +
  stat_compare_means(aes(group = Gender), label = "p.format", method = "kruskal.test")

filter: removed one row (<1%), 1,091 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.Med.Statin.LLD_byGender.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ LDLGroup, group.by = "Med.Statin.LLD", data = AERNASE.clin %>% filter(!is.na(LDLGroup) & !is.na(Med.Statin.LLD)), method = "kruskal.test")

filter: removed 458 rows (42%), 634 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(LDLGroup) & !is.na(Med.Statin.LLD)), 
                  x = c("LDLGroup"),
                  y = "HMOX1", 
                  xlab = "LDL (mg/dL) per LLD use",
                  ylab = "HMOX1 (normalized expression))",
                  color = "Med.Statin.LLD",
                  palette = c("#49A01D", "#1290D9"),
                  add = "jitter") +
  stat_compare_means(aes(group = Med.Statin.LLD), label = "p.format", method = "kruskal.test")

filter: removed 458 rows (42%), 634 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.LDLGroups_byMed.Statin.LLD.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ Hypercholesterolemia, group.by = "Med.Statin.LLD", data = AERNASE.clin %>% filter(!is.na(Hypercholesterolemia) & !is.na(Med.Statin.LLD)), method = "kruskal.test")

filter: removed 90 rows (8%), 1,002 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(Hypercholesterolemia) & !is.na(Med.Statin.LLD)), 
                  x = c("Hypercholesterolemia"),
                  y = "HMOX1", 
                  xlab = "Diagnosed hypercholesterolemia per LLD use",
                  ylab = "HMOX1 (normalized expression))",
                  color = "Med.Statin.LLD",
                  palette = c("#49A01D", "#1290D9"),
                  add = "jitter") +
  stat_compare_means(aes(group = Med.Statin.LLD), label = "p.format", method = "kruskal.test")

filter: removed 90 rows (8%), 1,002 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.LDLGroups_byMed.Statin.LLD.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

Kidney function (eGFR)

We want to create figures of CONVOCALS_downstream levels stratified by kidney function.

• Box and Whisker plot for CONVOCALS_downstream plaque levels by chronic kidney disease (CKD) group (1, 2, 3, 4, 5)
• Box and Whisker plot for CONVOCALS_downstream plaque levels by eGFR (MDRD-based) group (90+, 60-89, 30-59, <30)

library(dplyr)

AERNASE.clin <- AERNASE.clin %>% mutate(eGFRGroup = factor(case_when(GFR_MDRD < 15 ~ "<15",
                                                             GFR_MDRD >= 15  & GFR_MDRD <= 29 ~ "15-29",
                                                             GFR_MDRD >= 30  & GFR_MDRD <= 59 ~ "30-59",
                                                             GFR_MDRD >= 60  & GFR_MDRD <= 89 ~ "60-89",
                                                             GFR_MDRD >= 90 ~ "90+")))

table(AERNASE.clin$eGFRGroup, AERNASE.clin$Gender)

       
        female male
  <15        1    0
  15-29      3   11
  30-59     77  165
  60-89    165  409
  90+       47  155

table(AERNASE.clin$eGFRGroup, AERNASE.clin$KDOQI)

       
        No data available/missing Normal kidney function CKD 2 (Mild) CKD 3 (Moderate) CKD 4 (Severe) CKD 5 (Failure)
  <15                           0                      0            0                0              0               1
  15-29                         0                      0            0                0             14               0
  30-59                         0                      0            0              242              0               0
  60-89                         0                      0          574                0              0               0
  90+                           0                    202            0                0              0               0

Now we can draw some graphs of plaque CONVOCALS_downstream levels per sex and kidney function group as median ± interquartile range.

HMOX1

# Global test

compare_means(HMOX1 ~ eGFRGroup, data = AERNASE.clin %>% filter(!is.na(eGFRGroup)), method = "kruskal.test")

filter: removed 59 rows (5%), 1,033 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(eGFRGroup)), 
                  x = c("eGFRGroup"),
                  y = "HMOX1", 
                  xlab = "eGFR (mL/min per 1.73 m2)",
                  ylab = "HMOX1 (normalized expression)",
                  color = "eGFRGroup",
                  palette = "npg",
                  add = "jitter") +
  stat_compare_means(method = "kruskal.test")

filter: removed 59 rows (5%), 1,033 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.EGFR.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ eGFRGroup, group.by = "Gender",  data = AERNASE.clin %>% filter(!is.na(eGFRGroup)), method = "kruskal.test")

filter: removed 59 rows (5%), 1,033 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(eGFRGroup)), 
                  x = c("eGFRGroup"),
                  y = "HMOX1", 
                  xlab = "eGFR (mL/min per 1.73 m2) per gender",
                  ylab = "HMOX1 (normalized expression)",
                  color = "Gender",
                  palette = c("#D5267B", "#1290D9"),
                  add = "jitter") +
  stat_compare_means(aes(group = Gender), label = "p.format", method = "kruskal.test")

filter: removed 59 rows (5%), 1,033 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.EGFR_byGender.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ KDOQI, data = AERNASE.clin %>% filter(!is.na(KDOQI)), method = "kruskal.test")

filter: removed 28 rows (3%), 1,064 rows remaining

p1 <- ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(KDOQI)), 
                  x = c("KDOQI"),
                  y = "HMOX1", 
                  xlab = "Kidney function (KDOQI)",
                  ylab = "HMOX1 (normalized expression)",
                  color = "KDOQI",
                  palette = "npg",
                  add = "jitter") +
  stat_compare_means(aes(group = KDOQI), label = "p.format", method = "kruskal.test")

filter: removed 28 rows (3%), 1,064 rows remaining

ggpar(p1 + rotate_x_text(45), legend = "right") 
rm(p1)
ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.KDOQI.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ KDOQI, group.by = "Gender",   data = AERNASE.clin %>% filter(!is.na(KDOQI)), method = "kruskal.test")

filter: removed 28 rows (3%), 1,064 rows remaining

p1 <- ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(KDOQI)), 
                  x = c("KDOQI"),
                  y = "HMOX1", 
                  xlab = "Kidney function (KDOQI) per gender",
                  ylab = "HMOX1 (normalized expression)",
                  color = "Gender",
                  palette = c("#D5267B", "#1290D9"),
                  add = "jitter") +
  stat_compare_means(aes(group = Gender), label = "p.format", method = "kruskal.test")

filter: removed 28 rows (3%), 1,064 rows remaining

ggpar(p1 + rotate_x_text(45), legend = "right") 
rm(p1)
ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.KDOQI_byGender.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ eGFRGroup,  data = AERNASE.clin %>% filter(!is.na(eGFRGroup) & !is.na(KDOQI)), method = "kruskal.test")

filter: removed 59 rows (5%), 1,033 rows remaining

p1 <- ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(eGFRGroup) & !is.na(KDOQI)), 
                  x = c("eGFRGroup"),
                  y = "HMOX1", 
                  xlab = "eGFR (mL/min per 1.73 m2) by KDOQI group",
                  ylab = "HMOX1 (normalized expression)",
                  color = "KDOQI",
                  palette = "npg",
                  add = "jitter") +
  stat_compare_means(method = "kruskal.test")

filter: removed 59 rows (5%), 1,033 rows remaining

ggpar(p1, legend = "right")
rm(p1)
ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.EGFR_KDOQI.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

BMI

We want to create figures of CONVOCALS_downstream levels stratified by BMI.

• Box and Whisker plot for CONVOCALS_downstream plaque levels by BMI WHO group (underweight, normal, overweight, obese)
• Box and Whisker plot for CONVOCALS_downstream plaque levels by BMI group (<18.5, 18.5-24.9, 25, 29.9, 30-24.9, 35+)

library(dplyr)

AERNASE.clin <- AERNASE.clin %>% mutate(BMIGroup = factor(case_when(BMI < 18.5 ~ "<18.5",
                                                     BMI >= 18.5  & BMI < 25 ~ "18.5-24",
                                                     BMI >= 25  & BMI < 30 ~ "25-29",
                                                     BMI >= 30  & BMI < 35 ~ "30-35",
                                                     BMI >= 35 ~ "35+"))) 

# require(labelled)
# AERNASE.clin$BMI_US <- as_factor(AERNASE.clin$BMI_US)
# AERNASE.clin$BMI_WHO <- as_factor(AERNASE.clin$BMI_WHO)
# table(AERNASE.clin$BMI_WHO, AERNASE.clin$BMI_US)

table(AERNASE.clin$BMIGroup, AERNASE.clin$Gender)

         
          female male
  <18.5        6    4
  18.5-24    118  273
  25-29      120  373
  30-35       37   89
  35+         13   17

table(AERNASE.clin$BMIGroup, AERNASE.clin$BMI_WHO)

         
          No data available/missing Underweight Normal Overweight Obese
  <18.5                           0           9      0          0     0
  18.5-24                         0           0    391          0     0
  25-29                           0           0      0        492     0
  30-35                           0           0      0          0   126
  35+                             0           0      0          0    30

Now we can draw some graphs of plaque CONVOCALS_downstream levels per sex and age group as median ± interquartile range.

HMOX1

# Global test
compare_means(HMOX1 ~ BMIGroup,  data = AERNASE.clin %>% filter(!is.na(BMIGroup)), method = "kruskal.test")

filter: removed 42 rows (4%), 1,050 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(BMIGroup)), 
                  x = c("BMIGroup"),
                  y = "HMOX1", 
                  xlab = "BMI groups (kg/m2)",
                  ylab = "HMOX1 (normalized expression)",
                  # color = "Gender",
                  # palette = c("#D5267B", "#1290D9"),
                  color = "BMIGroup",
                  palette = "npg",
                  add = "jitter") +
  stat_compare_means(label = "p.format", method = "kruskal.test")

filter: removed 42 rows (4%), 1,050 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.BMI.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ BMIGroup, group.by = "Gender", data = AERNASE.clin %>% filter(!is.na(BMIGroup)), method = "kruskal.test")

filter: removed 42 rows (4%), 1,050 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(BMIGroup)), 
                  x = c("BMIGroup"),
                  y = "HMOX1", 
                  xlab = "BMI groups (kg/m2) per gender",
                  ylab = "HMOX1 (normalized expression)",
                  color = "Gender",
                  palette = c("#D5267B", "#1290D9"),
                  add = "jitter") +
  stat_compare_means(aes(group = Gender), label = "p.format", method = "kruskal.test")

filter: removed 42 rows (4%), 1,050 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.BMI_byGender.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ BMIGroup,  data = AERNASE.clin %>% filter(!is.na(BMIGroup) & !is.na(BMI_WHO)), method = "kruskal.test")

filter: removed 44 rows (4%), 1,048 rows remaining

p1 <- ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(BMIGroup) & !is.na(BMI_WHO)), 
                  x = c("BMIGroup"),
                  y = "HMOX1", 
                  xlab = "BMI groups (kg/m2) per WHO categories",
                  ylab = "HMOX1 (normalized expression)",
                  color = "BMI_WHO",
                  palette = "npg",
                  add = "jitter") +
  stat_compare_means(method = "kruskal.test")

filter: removed 44 rows (4%), 1,048 rows remaining

ggpar(p1, legend = "right")
rm(p1)
ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.BMI_byWHO.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

Diabetes

We want to create figures of CONVOCALS_downstream levels stratified by type 2 diabetes.

• Box and Whisker plot for CONVOCALS_downstream plaque levels by type 2 diabetes group (no, yes)

Now we can draw some graphs of plaque CONVOCALS_downstream levels per sex and age group as median ± interquartile range.

HMOX1

compare_means(HMOX1 ~ DiabetesStatus,  
              data = AERNASE.clin %>% filter(!is.na(DiabetesStatus)), method = "kruskal.test")

filter: no rows removed

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(DiabetesStatus)),
                  x = c("DiabetesStatus"),
                  y = "HMOX1",
                  xlab = "Diabetes status",
                  ylab = "HMOX1 (normalized expression)",
                  # color = "Gender",
                  # palette = c("#D5267B", "#1290D9"),
                  color = "DiabetesStatus",
                  palette = "npg",
                  add = c("median_iqr", "jitter")) +
  stat_compare_means(label = "p.format", method = "kruskal.test")

filter: no rows removed

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.Diabetes.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ DiabetesStatus, group.by = "Gender", data = AERNASE.clin %>% filter(!is.na(DiabetesStatus)), method = "kruskal.test")

filter: no rows removed

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(DiabetesStatus)),
                  x = c("DiabetesStatus"),
                  y = "HMOX1",
                  xlab = "Diabetes status per gender",
                  ylab = "HMOX1 (normalized expression)",
                  color = "Gender",
                  palette = c("#D5267B", "#1290D9"),
                  add = c("median_iqr", "jitter")) +
  stat_compare_means(aes(group = Gender), label = "p.format", method = "kruskal.test")

filter: no rows removed

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.Diabetes_byGender.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

NA
NA

Smoking

We want to create figures of CONVOCALS_downstream levels stratified by smoking.

• Box and Whisker plot for CONVOCALS_downstream plaque levels by smoking group (never, ex, current)

Now we can draw some graphs of plaque CONVOCALS_downstream levels per sex and age group as median ± interquartile range.

HMOX1

# Global test
compare_means(HMOX1 ~ SmokerStatus,  data = AERNASE.clin %>% filter(!is.na(SmokerStatus)), method = "kruskal.test")

filter: removed 45 rows (4%), 1,047 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(SmokerStatus)), 
                  x = c("SmokerStatus"),
                  y = "HMOX1", 
                  xlab = "Smoker status",
                  ylab = "HMOX1 (normalized expression)",
                  # color = "Gender",
                  # palette = c("#D5267B", "#1290D9"),
                  color = "SmokerStatus",
                  palette = "npg",
                  add = c("median_iqr", "jitter")) +
  stat_compare_means(label = "p.format", method = "kruskal.test")

filter: removed 45 rows (4%), 1,047 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.Smoking.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ SmokerStatus, group.by ="Gender", data = AERNASE.clin %>% filter(!is.na(SmokerStatus)), method = "kruskal.test")

filter: removed 45 rows (4%), 1,047 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(SmokerStatus)), 
                  x = c("SmokerStatus"),
                  y = "HMOX1", 
                  xlab = "Smoker status per gender",
                  ylab = "HMOX1 (normalized expression)",
                  color = "Gender",
                  palette = c("#D5267B", "#1290D9"),
                  add = c("median_iqr", "jitter")) +
  stat_compare_means(aes(group = Gender), label = "p.format", method = "kruskal.test")

filter: removed 45 rows (4%), 1,047 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.Smoking_byGender.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

Stenosis

We want to create figures of CONVOCALS_downstream levels stratified by stenosis grade.

• Box and Whisker plot for CONVOCALS_downstream plaque levels by stenosis grade group (<70, 70-89, 90+)

library(dplyr)

AERNASE.clin <- AERNASE.clin %>% mutate(StenoticGroup = factor(case_when(stenose == "0-49%" ~ "<70",
                                                     stenose == "0-49%" ~ "<70",
                                                     stenose == "50-70%" ~ "<70",
                                                     stenose == "70-90%" ~ "70-89",
                                                     stenose == "50-99%" ~ "90+",
                                                     stenose == "70-99%" ~ "90+",
                                                     stenose == "100% (Occlusion)" ~ "90+",
                                                     stenose == "90-99%" ~ "90+")))

table(AERNASE.clin$StenoticGroup, AERNASE.clin$Gender)

       
        female male
  <70       17   61
  70-89    151  358
  90+      133  342

table(AERNASE.clin$stenose, AERNASE.clin$StenoticGroup)

                  
                   <70 70-89 90+
  missing            0     0   0
  0-49%              5     0   0
  50-70%            73     0   0
  70-90%             0   509   0
  90-99%             0     0 438
  100% (Occlusion)   0     0  10
  NA                 0     0   0
  50-99%             0     0   4
  70-99%             0     0  23
  99                 0     0   0

Now we can draw some graphs of plaque CONVOCALS_downstream levels per sex and age group as median ± interquartile range.

HMOX1

# Global test
compare_means(HMOX1 ~ StenoticGroup,  data = AERNASE.clin %>% filter(!is.na(StenoticGroup)), method = "kruskal.test")

filter: removed 30 rows (3%), 1,062 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(StenoticGroup)), 
                  x = c("StenoticGroup"),
                  y = "HMOX1", 
                  xlab = "Stenotic grade",
                  ylab = "HMOX1 (normalized expression)",
                  # color = "Gender",
                  # palette = c("#D5267B", "#1290D9"),
                  color = "StenoticGroup",
                  palette = "npg",
                  add = "jitter") +
  stat_compare_means(label = "p.format", method = "kruskal.test")

filter: removed 30 rows (3%), 1,062 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.Stenosis.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

compare_means(HMOX1 ~ StenoticGroup, group.by = "Gender", data = AERNASE.clin %>% filter(!is.na(StenoticGroup)), method = "kruskal.test")

filter: removed 30 rows (3%), 1,062 rows remaining

ggpubr::ggboxplot(AERNASE.clin %>% filter(!is.na(StenoticGroup)), 
                  x = c("StenoticGroup"),
                  y = "HMOX1", 
                  xlab = "Stenotic grade per gender",
                  ylab = "HMOX1 (normalized expression)",
                  color = "Gender",
                  palette = c("#D5267B", "#1290D9"),
                  add = "jitter") +
  stat_compare_means(aes(group = Gender), label = "p.format", method = "kruskal.test")

filter: removed 30 rows (3%), 1,062 rows remaining

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.Stenosis_byGender.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

Symptoms

We want to create per-symptom figures.

library(dplyr)

table(AERNASE.clin$AgeGroup, AERNASE.clin$AsymptSympt2G)

       
        Asymptomatic Symptomatic
  <55             10          44
  55-64           32         245
  65-74           47         385
  75-84           20         277
  85+              2          30

table(AERNASE.clin$Gender, AERNASE.clin$AsymptSympt2G)

        
         Asymptomatic Symptomatic
  female           24         284
  male             87         697

table(AERNASE.clin$AsymptSympt2G)

Asymptomatic  Symptomatic 
         111          981 

Now we can draw some graphs of plaque CONVOCALS_downstream levels per symptom group as median ± interquartile range.

HMOX1

# ?ggpubr::ggboxplot()
my_comparisons <- list(c("Asymptomatic", "Symptomatic"))

p1 <- ggpubr::ggboxplot(AERNASE.clin, 
                  x = "AsymptSympt2G", y = "HMOX1",
                  title = "HMOX1 (normalized expression) levels per symptom", 
                  xlab = "Symptoms",
                  ylab = "HMOX1 (normalized expression)",
                  color = "AsymptSympt2G", 
                  # palette = c(uithof_color[16], uithof_color[23]),
                  palette = "npg",
                  add = "dotplot", # Add dotplot
                  add.params = list(binwidth = 0.1, dotsize = 0.3)
          ) +
  stat_compare_means(comparisons = my_comparisons, method = "wilcox.test")
ggpar(p1, legend = c("right"), legend.title = "Symptoms")

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.AsymptSympt2G.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

rm(p1)

compare_means(HMOX1 ~ AsymptSympt2G, group.by = "Gender", data = AERNASE.clin, method = "kruskal.test")

p1 <- ggpubr::ggboxplot(AERNASE.clin, 
                  x = "AsymptSympt2G", y = "HMOX1",
                  title = "HMOX1 (normalized expression) levels per symptom by gender", 
                  xlab = "Symptoms",
                  ylab = "HMOX1 (normalized expression)",
                  color = "Gender",
                  palette = c("#D5267B", "#1290D9"),
                  add = "dotplot", # Add dotplot
                  add.params = list(binwidth = 0.1, dotsize = 0.3)
          ) +
  stat_compare_means(aes(group = Gender), label = "p.format",  method = "wilcox.test")
ggpar(p1, legend = c("right"), legend.title = "Symptoms")

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.AsymptSympt2G.byGender.pdf"), plot = last_plot())

Saving 7.29 x 4.51 in image

rm(p1)

Forest plots

We would also like to visualize the multivariable analyses results.

library(ggplot2)
library(openxlsx)
model1_target <- read.xlsx(paste0(OUT_loc, "/", Today, ".AERNASE.clin.targets.Bin.Uni.",TRAIT_OF_INTEREST,".RANK.Symptoms.MODEL1.xlsx"))
model2_target <- read.xlsx(paste0(OUT_loc, "/", Today, ".AERNASE.clin.targets.Bin.Multi.",TRAIT_OF_INTEREST,".RANK.Symptoms.MODEL2.xlsx"))
model1_target$model <- "univariate"
model2_target$model <- "multivariate"

models_target <- rbind(model1_target, model2_target)
models_target

NA

Forest plots.

HMOX1

dat <- data.frame(group = factor(c("Age, sex-adjusted", "Age, sex, and adjusted for risk factors"), 
                           
                           levels=c("Age, sex, and adjusted for risk factors", "Age, sex-adjusted")),
                  cen = c(models_target$OR[models_target$Predictor=="HMOX1"]),
                  low = c(models_target$low95CI[models_target$Predictor=="HMOX1"]),
                  high = c(models_target$up95CI[models_target$Predictor=="HMOX1"]))

fp <- ggplot(data = dat, aes(x = group, y = cen, ymin = low, ymax = high)) +
  geom_pointrange(linetype = 2, size = 1, colour = c("#1290D9", "#49A01D")) + 
  geom_hline(yintercept = 1, lty = 2) +  # add a dotted line at x=1 after flip
  coord_flip(ylim = c(0.8, 1.7)) +  # flip coordinates (puts labels on y axis)
  xlab("Model") + ylab("OR (95% CI) for symptomatic plaques") +
  ggtitle("Plaque HMOX1 normalized expression (1 SD increment, n = 622)") +
  theme_minimal()  # use a white background
print(fp)

ggsave(file = paste0(PLOT_loc, "/",Today,".AERNASE.clin.",TRAIT_OF_INTEREST,".HMOX1.plaque.forest.pdf"), plot = fp)

Saving 7.29 x 4.51 in image

rm(fp)

Target expression vs. cytokines plaque levels correlations

We will plot the correlations of other cytokine plaque levels to the CONVOCALS_downstream plaque levels. These include:

• IL2
• IL4
• IL5
• IL6
• IL8
• IL9
• IL10
• IL12
• IL13
• IL21
• INFG
• TNFA
• MIF
• MCP1
• MIP1a
• RANTES
• MIG
• IP10
• Eotaxin1
• TARC
• PARC
• MDC
• OPG
• sICAM1
• VEGFA
• TGFB

In addition we will look at three metalloproteinases which were measured using an activity assay.

• MMP2
• MMP8
• MMP9

The proteins were measured using FACS and LUMINEX. Given the different platforms used (FACS vs. LUMINEX), we will inverse rank-normalize these variables as well to scale them to the same scale as the CONVOCALS_downstream` plaque levels.

We will set the measurements that yielded ‘0’ to NA, as it is unlikely that any protein ever has exactly 0 copies. The ‘0’ yielded during the experiment are due to the limits of the detection.

Prepare data

# fix names
names(AEDB.CEA)[names(AEDB.CEA) == "VEFGA"] <- "VEGFA"

# fix names
names(AERNASE.clin)[names(AERNASE.clin) == "IL6"] <- "IL6rna"
names(AERNASE.clin)[names(AERNASE.clin) == "MMP9"] <- "MMP9rna"

cytokines <- c("IL2", "IL4", "IL5", "IL6", "IL8", "IL9", "IL10", "IL12", "IL13", "IL21", 
               "INFG", "TNFA", "MIF", "MCP1", "MIP1a", "RANTES", "MIG", "IP10", "Eotaxin1", 
               "TARC", "PARC", "MDC", "OPG", "sICAM1", "VEGFA", "TGFB")
metalloproteinases <- c("MMP2", "MMP8", "MMP9")


AERNASE.clin <- merge(AERNASE.clin, 
                            subset(AEDB.CEA, select = c("STUDY_NUMBER", 
                                                        cytokines,
                                                        metalloproteinases)), 
                            by.x = "study_number", by.y = "STUDY_NUMBER", sort = TRUE, all.x = TRUE)

proteins_of_interest <- c(cytokines, metalloproteinases)

proteins_of_interest_rank = unlist(lapply(proteins_of_interest, paste0, "_rank"))

# make variables numerics()
AERNASE.clin <- AERNASE.clin %>%
  mutate_each(funs(as.numeric), proteins_of_interest)

Warning: `funs()` was deprecated in dplyr 0.8.0.
Please use a list of either functions or lambdas: 

  # Simple named list: 
  list(mean = mean, median = median)

  # Auto named with `tibble::lst()`: 
  tibble::lst(mean, median)

  # Using lambdas
  list(~ mean(., trim = .2), ~ median(., na.rm = TRUE))Warning: `mutate_each()` was deprecated in dplyr 0.7.0.
Please use `across()` instead.

Error in `mutate_each_impl()`:
! Can't select columns that don't exist.
✖ Columns `IL2`, `IL4`, `IL5`, `IL8`, `IL9`, etc. don't exist.
Backtrace:
 1. AERNASE.clin %>% ...
 2. dplyr::mutate_each(., funs(as.numeric), proteins_of_interest)
 3. dplyr:::mutate_each_impl(tbl, funs, enquos(...), "mutate_each")

Visualize transformations

We will just visualize these transformations.

Correlations

Here we calculate correlations between CONVOCALS_downstream and 28 other cytokines. We use Spearman’s test, thus, correlations a given in rho. Please note the indications of measurement methods:

• L: LUMINEX
• E: ELISA
• a: activity assay

While visually attractive we are not necessarily interested in the correlations between all the cytokines, rather of CONVOCALS_downstream` with other cytokines only.

HMOX1

Another version - probably not good.

Target expression vs. cytokines plaque levels lm()

Model 1

In this model we correct for Age, Gender, and year of surgery.

Here we use the inverse-rank normalized data - visually this is more normally distributed.

Analysis of plaque cytokines traits as a function of plaque CONVOCALS_downstream levels.

Model 2

In this model we correct for Age, Gender, year of surgery, Hypertension status, Diabetes status, current smoker status, lipid-lowering drugs (LLDs), antiplatelet medication, eGFR (MDRD), BMI, MedHx_CVD (combination of CAD history, stroke history, and peripheral interventions), and stenosis.

Here we use the inverse-rank normalized data - visually this is more normally distributed.

Analysis of plaque cytokines as a function of plaque CONVOCALS_downstream levels.

DT::datatable(GLM.results)

# Save the data
cat("Writing results to Excel-file...\n")
### Univariate
library(openxlsx)
write.xlsx(GLM.results,
           file = paste0(OUT_loc, "/",Today,".AERNASE.clin.Con.Multi.",TRAIT_OF_INTEREST,"_Plaque.Cytokines_Plaques.RANK.MODEL2.xlsx"),
           rowNames = FALSE, colNames = TRUE, sheetName = "Con.Multi.PlaquePheno")
# Removing intermediates
cat("Removing intermediate files...\n")
rm(TRAIT, trait, currentDF, GLM.results, GLM.results.TEMP, fit, model_step)

Target expression vs. vulnerability index

Here we plot the levels of inverse-rank normal transformed CONVOCALS_downstream plaque levels from experiment 1 and 2 to the Plaque vulnerability index.

Visualisations

HMOX1

Model 1

In this model we correct for Age, Gender, and year of surgery.

Here we use the inverse-rank normalized data - visually this is more normally distributed.

Analysis of the plaque vulnerability index as a function of plaque CONVOCALS_downstream levels.

Model 2

In this model we correct for Age, Gender, Hypertension status, Diabetes status, current smoker status, lipid-lowering drugs (LLDs), antiplatelet medication, eGFR (MDRD), BMI, MedHx_CVD (combination of CAD history, stroke history, and peripheral interventions), and stenosis..

Session information

---

Version:      v1.1.1
Last update:  2024-01-09
Written by:   Sander W. van der Laan (s.w.vanderlaan-2[at]umcutrecht.nl).
Description:  Script to analyse Targets from the Ather-Express Biobank Study.
Minimum requirements: R version 3.5.2 (2018-12-20) -- 'Eggshell Igloo', macOS Mojave (10.14.2).

**MoSCoW To-Do List**
The things we Must, Should, Could, and Would have given the time we have.
_M_

_S_

_C_

_W_

**Changes log**
* v1.1.1 Textual fixes.
* v1.1.0 Update to study database; update to bulk RNAseq data (deeper sequenced).
* v1.0.1 Fix to the start of this notebook.
* v1.0.0 Inital version.

---

Saving environment

save.image(paste0(PROJECT_loc, "/",Today,".",PROJECTNAME,".bulkRNAseq.additional_figures.RData"))

© 1979-2024 Sander W. van der Laan | s.w.vanderlaan[at]gmail[dot]com | vanderlaanand.science.